/* Table 1 */
/* Basic models, with timing and specification tests */
/* be sure mixlogit package is installed */

#delimit;
set mem 30m;
set more off; 

log using Table1.log, t replace;

use formation_new;

timer on 1;

clogit realg anmax2 minor minwin numpar median mginvest pro anti dompar sq prevpm mgodiv1 gdiv1, group(case); 

timer off 1;
timer list 1;

/* Lagrangian test */

predict cp, p;

sort case;

by case: egen anmax2S = sum(anmax2 * cp);
by case: egen minorS = sum(minor * cp);
by case: egen minwinS = sum(minwin * cp);
by case: egen numparS = sum(numpar * cp);
by case: egen medianS = sum(median * cp);
by case: egen mginvestS = sum(mginvest * cp);
by case: egen proS = sum(pro * cp);
by case: egen antiS = sum(anti * cp);
by case: egen domparS = sum(dompar * cp);
by case: egen sqS = sum(sq * cp);
by case: egen prevpmS = sum(prevpm * cp);
by case: egen mgodiv1S = sum(mgodiv1 * cp);
by case: egen gdiv1S = sum(gdiv1 * cp);

generate anmax2Z = 0.5*(anmax2 - anmax2S)^2;
generate minorZ = 0.5*(minor - minorS)^2;
generate minwinZ = 0.5*(minwin - minwinS)^2;
generate numparZ = 0.5*(numpar - numparS)^2;
generate medianZ = 0.5*(median - medianS)^2;
generate mginvestZ = 0.5*(mginvest - mginvestS)^2;
generate proZ = 0.5*(pro - proS)^2;
generate antiZ = 0.5*(anti - antiS)^2;
generate domparZ = 0.5*(dompar - domparS)^2;
generate sqZ = 0.5*(sq - sqS)^2;
generate prevpmZ = 0.5*(prevpm - prevpmS)^2;
generate mgodiv1Z = 0.5*(mgodiv1 - mgodiv1S)^2;
generate gdiv1Z = 0.5*(gdiv1 - gdiv1S)^2;

/* McFadden-Train LM test */

clogit realg 
anmax2 minor minwin numpar median mginvest pro anti dompar sq prevpm mgodiv1 gdiv1
anmax2Z minorZ minwinZ numparZ medianZ mginvestZ proZ antiZ domparZ sqZ prevpmZ mgodiv1Z gdiv1Z, 
group(case); 

/* significant: anmax2 numpar median dompar sq.  t>1: minor, mgodiv1 */
/* anmax2 barely identified even as fixed coefficient, so don't include as random */

replace gdiv1 = gdiv1 * -1; 

/* all random (except three variables w/ separation problem) */

timer on 2;

mixlogit realg anmax2 pro anti, rand(minor minwin mginvest prevpm numpar dompar sq median mgodiv1 gdiv1) ln(3) nrep(200) group(case);

timer off 2;
timer list 2;

/* generate random start values */

matrix coeffs = e(b);
matrix stderrs = vecdiag(e(V));
matrix direct = matuniform(1,23) - J(1,23,0.5);
matrix b0 = coeffs + direct;

local i = 1;

while `i' < 6 {;

  matrix list b0;

  mixlogit realg anmax2 pro anti, rand(minor minwin mginvest prevpm numpar dompar sq median mgodiv1 gdiv1) ln(3) nrep(200) group(case) from(b0, copy);

  matrix direct = (matuniform(1,23) - J(1,23,0.5)) * `i';
  matrix b0 = coeffs + direct;

  local i = `i' + 1;

};


/* specification based on Lagrange multiplier test */

replace gdiv1 = gdiv1*-1;

timer on 3;

mixlogit realg anmax2 pro anti minwin mginvest prevpm gdiv1, rand(minor numpar dompar sq median mgodiv1) ln(2) nrep(5) group(case); 

timer off 3;
timer list 3;

/* generate random start values */

matrix coeffs = e(b);
matrix stderrs = vecdiag(e(V));
matrix direct = matuniform(1,19) - J(1,19,0.5);
matrix b0 = coeffs + direct;

local k = 1;

while `k' < 6 {;

  matrix list b0;

  mixlogit realg anmax2 pro anti minwin mginvest prevpm gdiv1, rand(minor numpar dompar sq median mgodiv1) ln(2) nrep(200) group(case)  from(b0, copy);

  matrix direct = (matuniform(1,19) - J(1,19,0.5)) * `k';
  matrix b0 = coeffs + direct;

  local k = `k' + 1;

};
